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ABSTRACT 



Context. Solar active regions are formed through the emergence of magnetic flux from the deeper convection zone. Recent satellite 
observations have shown that a horizontal divergent flow (HDF) stretches out over the solar surface just before the magnetic flux 
appearance. 

Aims. The aims of this study are to investigate the driver of the HDF and to see the dependency of the HDF on the parameters of the 
magnetic flux in the convection zone. 

Methods. We conduct three-dimensional magnetohydrodynamic (3D MHD) numerical simulations of the magnetic flux emergence 
and vary the parameters in the initial conditions. An analytical approach is also taken to explain the dependency. 
Results. The horizontal gas pressure gradient is found to be the main driver of the HDF. The maximum HDF speed shows positive 
correlations with the field strength and twist intensity. The HDF duration has a weak relation with the twist, while it shows negative 
dependency on the field strength only in the case of the stronger field regime. 

Conclusions. Parametric dependencies analyzed in this study may allow us to probe the structure of the subsurface magnetic flux by 
observing properties of the HDF. 



1. Introduction 

The dynamics of the rising magnetic flux and the forma- 
tion process of solar active regions have widely been in- 
vestigated through a series of analytical and numerical stud- 
ies. Parked (Il975l) first calculated the rising speed of a flux 
tube in the convection zone (CZ) by considering a force bal- 
ance between the magnetic buoyanc y and the aero dynamic 
drag acting on the flux tube, while Schiissler (1979) simu- 
lated a buoyant emergence of a flux tube in the CZ in a two- 
dim ensional (2D) scheme. iM oreno-Ins ertis & Emonel fl996) 
and lEmonet & Moreno-Insertisl (1 19981) investigated the depen- 
dency of the rising tube on the initial twist intensity and found 
that th e tube needs a certa in degree of twis t to ho ld its coherency. 
After IShibata et al.l (1 19891) and iMagaral (1200 ll) conducted 2D 
simulations of flux emergence from the surface layer into the 
cor ona, 3D sim ula tions have been carrie d out by, among oth- 
ers. lFa3 d200lh and lArchontis et al.l d2004l) . lMurrav etaD(l2006h 
simulated other 3D flux emergences of this kind and surveyed 
the tube's dependency on the initial field strength and twist in- 
tensity. 

Recently, iToriumi & Yokovamal (I201QL l201ll 120121) com- 
bined the CZ, the photosphere/chromosphere, and the corona 
into a single computational domain and simulated magnetic flux 
emergence from a deeper CZ both in 2D and 3D. As a result, the 
initial flux placed at a depth of -20 Mm starts its emergence in 
the solar interior, which then slows down gradually in the upper- 
most CZ. This is because the plasma pushed up by the emerging 
flux rises to the isothermally- stratified (i.e. convectively- stable) 
surface layer and is then trapped and compressed between them, 



which, in turn, suppresses the rising flux from below. Such com- 
pressed plasma will escape laterally around the photospheric 
layer from the rising flux as a horizontal divergent flow (HDF), 
just before the flux itself reaches the surface. Using SDO/HMI 
data, IToriumi et al.l (1201 2l) observed the emerging active region 
located away from the solar disk center and found the HDF in 
the Dopplergram, up to about 100 min before the start of the 
flux emergence. 

In the present study, we report the results of the parametric 
survey of the 3D magnetohydrodynamic (MHD) flux emergence 
simulation. The aims of this study are to investigate which force 
drives the HDF and to observe the dependence of the HDF on 
the parameters in the simulation. One important feature of this 
HDF study is that it can be a probe for exploring the physical 
state of the magnetic field in the upper CZ. That is, we may be 
able to obtain valuable information on the subsurface layers from 
the direct optical observation at the surface. Therefore, in this 
numerical study, we vary the parameters of the initial flux tube, 
and then check the characteristics of the consequent HDF seen 
at the surface layer. 

In the next section we introduce the basic setup of the numer- 
ical calculation. In Section[3]we show the results of the paramet- 
ric survey, and in Section [4] we provide some analytic explana- 
tions of the results. We finally summarize the paper in Section 




2. Numerical Setup 

The basic MHD equations, normalizing units, computational do- 
main size, grid spacings, boundary conditions, and background 
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strati fication are the same as those in iToriumi & Yokoyamal 
(2012). The MHD equations in vector form are: 



dp 
dt 



+ V • (pV) = 0, 



|(pV) + V.(pW + p/-£ + £/|-p* = 0, 



dB 

dt 



V x (V x B), 



(1) 



(2) 



(3) 



+V- 

and 
U = 



Lu + p+±pV 2 \v+±ExB 



Pg-V = 0, 



1 



y-\p 



1 



VxB, 



(4) 



(5) 



(6) 



fa „ 

P = —pT, 
m 



where p denotes the gas density, V velocity vector, p pressure, 
B magnetic field, c the speed of light, E electric field, and T 
temperature, while U is the internal energy per unit mass, / the 
unit tensor, fa the Boltzmann constant, m(= const.) the mean 
molecular mass, and g the uniform gravitational acceleration. 
We assume the medium to be an inviscid perfect gas with a spe- 
cific heat ratio y = 5/3. All the physical values are normalized 
by the pressure scale height Ho = 200 km for length, the sound 
speed C s o = 8 km s" 1 for velocity, To = Ho/C s o - 25 s for 
time, and po = 1.4 x 10" 7 g cm -3 for density, all of which are the 
typical values in the photosphere. The units for pressure, temper- 
ature, and magnetic field strength are po = 9.0 x 10 4 dyn cm -2 , 
T = 4000 K, and B = 300 G, respectively. 

Here, 3D Cartesian coordinates (x, y, z) are used, where z is 
parallel to the gravitational acceleration vector, g = (0, 0, -go), 
and go = C^JiyHo) by definition. The simulation domain is 
(-400, -200, -200) < (x/H , y/H , z/H ) < (400, 200, 250), re- 
solved by 1602 x 256 x 1024 grids. In the x-direction, the mesh 
size is Ax/ Ho = 0.5 (uniform). In the y-direction (z-direction), 
the mesh size is Ay /Ho = 0.5 {Az/Hq = 0.2) in the central area 
of the domain, which gradually increases for each direction. We 
assume periodic boundaries for both horizontal directions and 
symmetric boundaries for the vertical direction. 

The background atmosphere consists of three different lay- 
ers. From the bottom, the layers are the adiabatically stratified 
CZ, the cool isothermal photo sphere/chromo sphere, and the hot 
isothermal corona. The stratification in the CZ (z/Hq < 0) is 
given as 



T = T ph -z 



dT 



dz 



(8) 



where T^/To - 1 is the respective temperature in the photo- 
sphere/chromosphere and 



dT 



dz 



y ~ 1 mgo 

y fa 



(9) 



Table 1. Summary of the simulation cases 



Case fl 


Field strength B tuhe /B 


Twist qH 


Wavelength A/H 


A 


67 


0.1 


400 


B 


133 


0.1 


400 


C 


33 


0.1 


400 


D 


67 


0.2 


400 


E 


67 


0.05 


400 


F 


67 


0.1 


100 


G 


67 


0.1 


25 



Note s, {a) Case A is the same as that simulated in lToriumi & Y okoyama 
(120121) . Cases B and C are for different field strengths than that of case 
A, while D and E are for different twists, and F and G different wave- 
lengths. 

is the adiabatic temperature gradient. The profile above the sur- 
face is 



jtanh 


Z Zcor 






W tr 





(10) 



where T COY /To = 100 is the temperatures in the corona, z C or/#o = 
10 is the base of the corona, and w tr /Ho = 0.5 is the transition 
scale length. Based on the temperature distribution above, the 
pressure and density profiles are defined by the equation of static 
pressure balance: 

dp(z) 



(J) dz 



+ p(z)go = 0. 



(11) 



The initial flux tube is embedded in the CZ at Ztube/^o = 
-100, i.e., Ztube = -20 Mm, of which the axial and azimuthal 
profiles are given as 

r 2 

B x (r) = £ tube exp| 



B+(f) = qrB x {r) 



p2 
^tube 



(12) 



respectively, where # t ube is the axial field strength, r the radial 
distance from the tube's center Cy tu be/#o,£tube/#o) = (0, -100), 
J?tube me typical radial size, and q the twist intensity. For the 
pressure balance between the field and the plasma, the pressure 
distribution inside the tube is defined as p x - p(z) + Sp QXC (the 
subscript "i" denotes inside the tube), where the pressure excess 
dPcxc(< 0) is described as 

B * (r) U^-A-i}. as) 



8tt 



p(z) + Sp e 



The density inside the tube is also defined as p x 
where 

^c=P(z)^exp(4), (14) 

and A is the perturbation wavelength. That is, the middle of the 
tube, x/Ho = 0, is in thermal equilibrium with external me- 
dia and is most buoyant. The buoyancy decreases as \x\/Hq in- 
creases. 

The parameters we varied are the field strength # t ube, 
the twist q, and the perturbation wavelength X. Table [T] 
summarizes the case s in th is study. The case simulated in 
IToriumi & Yokoyamal (120121) is named here as case A, while 
cases B-F are for different field strength, twist, and wavelength 
than those of A. Here we fixed the tube's radial size at R^q/Ho = 
5 for all the cases. It should be noted that the crit ical twist for the 
kink instability is qHo = 0.2 dLinton et al.ll 19961) . Therefore, all 
the tubes examined here are stable or, at least, marginally stable 
against the instability at the beginning of the calculation. 
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Fig. 1. Height- time evolution of the flux tube, (a) Cases for dif- 
ferent #tube- (b) Cases for different q. (c) Cases for different A. 



3. Simulation Results 

3.1. General Evolution 

Fig. Q] shows the temporal evolution of the apex of the rising 
tube, Ztop(0- Also, in Fig. [2 we plot the total field strength, 
log 10 (\B\/Bq), of the initial condition for case A and the final 
states for all the cases. As can be seen in Fig. [T^, it is clear 
that the tubes with stronger field # t ube rise faster. The rising 
speed of each tube in the CZ is in simple proporti on to the ini- 
tial fie ld strength, which is well in accordance with Mu rray et al.l 
(2006) and other previous studies. Case A, which has a middle 
field strength, shows deceleration just before it reaches the sur- 
face. This deceleration is the result of the plasma accumulation, 
which is caused by the trapping of material between the rising 
tube and the isothermally- stratified photosphere above the tube. 
After a while, the tube then starts further emergence into the at- 
mosphere (see the top panels in Fig.0. As for the strongest case 
in B, the accumulation becomes less marked, and thus the tube 
almost directly passes through the surface layer and expands into 
the higher corona, without undergoing strong deceleration (Fig. 
03). When the field is very weak, as in case C, the tube stops its 
emergence halfway to the surface, since the tube's buoyancy is 
not strong enough to continue its emergence (Fig.[2t). 

Fig. [TJ) shows the evolution with different twist q. In this 
figure, all three tubes are seen to rise almost at the same rate 
in the CZ, which is again consistent with previous studies (e.g. 




Iog 10 (l5l/5 ) 



Fig. 2. Total magnetic field strength of the initial condition for 
case A and the final states for cases A to G, plotted over the 
range -400 < x/H < and < y/H < 200. 



Murrav et al. 2006). When the twist is weak and thus cannot hold 
the coherency (case E), the tube expands and suffers deformation 
by aerodynamic drag. As a result, the tube cannot maintain a 
strong enough magnetic field to continue its further emergence 
(Fig. EE). 

Fig. [J: compares three cases with different wavelengths of 
the initial perturbation. The initial wavelength is crucial for 
two factors: the curvature force (magnetic tension), which pulls 
down the rising tube, and the drainage of the internal media due 
to gravity, which encourages the emergence. When the wave- 
length is smaller, the curvature force is expected to be stronger, 
while the drainage becomes more effective. In Fig.Q]: the rising 
velocities of cases A and F are almost the same, which indicates 
that both effects cancel each other out. However, the shortest 
wavelength tube (case G) shows a much slower emergence rate 
in the CZ, which indicates that the curvature force is more effec- 
tive and slows down the emergence. As for the emergence above 
the surface layer, on the contrary, Fig. OJ shows an exactly op- 
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Fig. 3. (a) Cross-sectional profile of the rising magnetic flux 
tube (case A). Plotted value is the logarithmic field strength 
log lo (\B\/Bo) averaged over 6.5 < x/Hq < 13.5, at the time 
t/ T Q = 500. The color saturates at \B\/B = 1.0 (red) and -4.0 
(purple), (b) Horizontal velocity V y /C s o (thick solid), pressure 
gradient -dp/dy x (lQHo/po) (dashed), and magnetic pressure 
gradient -dp m /dy x (10H /p ) (dash-dotted), at t/r = 500. (c) 
Same as (b) but for t/ro = 600. 



posite trend that the rising is much faster when the wavelength 
is shortest (case G). One may find that, in Fig[2fj, the main tube 
remains in the CZ at around z/Hq = -50, while the upper part 
has been detached from the main tube and has started further 
emergence into the atmosphere. The reason for the rapid ascent 
may be because, in the shortest wavelength case, namely, in the 
highly curved loops, the draining of the plasma from the apex 
is more effective, which helps the faster emergence above the 
photosphere. 

3.2. Driver of the HDF 

Fig. [3^ is the cross-sectional distribution of the field strength 
of Case A. The plotted value is the logarithmic field strength 
log 10 (\B\/Bq) averaged over 6.5 < x/Hq < 13.5, i.e., around 
x/Ho = 10. The reason we choose this x-range is to se- 
lect one folded struct ure at the tube's surface (see Fig. 3 of 
iToriumi & Yokoyamall2012l) . In this figure, there is a flow field 
in front of the rising flux tube, which is flowing from the apex 
to the flanks of the tube. One characteristic of this plasma layer 
is the horizontal divergent flow (HDF) that is seen at the solar 
surface just before the flux tube itself emerges. 



To investigate which force drives the HDF, in Figs. [3j) 
and c, we plot the horizontal flow velocity V y , pressure gra- 
dient -dp/dy, and magnetic pressure gradient -dp m /dy, aver- 
aged over 6.5 < x/H < 13.5 and -10 < z/H < 0, where 
p m = B 2 /(8n). Magnetic tension is not plotted here, since it is 
rather small compared to the two other forces. At t/ro = 500, 
before the tube reaches the uppermost CZ, z/Hq > -10, the 
horizontal flow is clearly driven only by the gas pressure, and, 
of course, the magnetic pressure gradient is zero. Therefore, we 
can conclude that the HDF prior to the flux appearance is caused 
by the pressure gradient. This is consistent with other numerica l 
simulations including thermal convection (iCheung et al.ll20ldl) . 
At t/ro = 600, the shallow layer is covered by the rising tube and 
the gas pressure gradient reverses its sign. Instead, the magnetic 
pressure gradient becomes dominant enough to drive the flow. 

3.3. Dependence of the HDF 

In this subsection, we show the dependence of the HDF on the 
initial field strength # tu b e and on the twist q. The investigated 
parameters are the duration of the HDF (from the HDF start 
to flux appearance), At/ To, and the maximum HDF velocity, 
max (Vy)/C S Q, during this time period. Here we defined the start 
time of the HDF as "when the horizontal speed V y in the horizon- 
tal range -50 < y/Ho < 50, averaged over 6.5 < x/Hq < 13.5 
and -10 < z/H < 0, exceeds 0.06C S (= 0.5 kms" 1 )" and the 
flux appearance as "when the field strength \B\ in this range ex- 
ceeds 0.67£ (= 200 G)." 

Figs. @k and b show the dependence of the duration on the 
field strength # tu b e and the twist q. Panel (a) is the comparison 
among the different field strength cases. A comparison of cases 
A and B, the middle and stronger field tubes, shows that the time 
duration is longer for the stronger field. If other stronger cases 
are considered (here we also plot two stronger tube cases other 
than A, B, and C), however, it may be found that the duration 
decreases with field strength. Thus, we can divide these cases 
into two groups: stronger cases that show a decreasing trend, 
which is fitted by a function of # t ~ be , and a middle case that 
deviates from the decreasing trend. The weakest tube, case C, 
did not reach the surface. That is why the duration is for case 
C. In contrast, in Panel (b), the duration is almost constant for 
the different twist cases. 

Dependence of the maximum HDF speed, max(V y ), is 
shown in Figs. ^ and d. Panel (c) indicates the positive linear 
correlation with the field strength # t ube- Again, the speed of case 
B is plotted as zero, since it did not reach the surface. Note that, 
in Panel (a), we found a gap between the middle-field regime 
and the stronger- field regime. Thus the linear fitting in Panel (c) 
might not reflect the actual trend. Nevertheless, the maximum 
speed basically increases with field strength. In Panel (d), we 
can see that the maximum HDF velocity is clearly proportional 
to the initial twist q. 



4. Analytic Explanation 

In this section, the dependencies of the rising speed and of the 
HDF on the physical parameters obtained in Section [3] are ana- 
lytically explained. 



4.1. Rising Speed 



In Section [3711 we found that the rising speed of the flux tube is 
proportional to the field strength # tu be an d the dependence on the 
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Fig. 4. (a) Dependence of the HDF duration 
on the field strength 5 tube , and (b) that on 
the twist q. (c) Dependence of of the maxi- 
mum HDF speed max (V y ) on the field strength 
Btube, and (d) that on the twist q. In Panel 
(a) and (c) we also plot other field strength 
cases, which are indicated by smaller asterisks. 
Solid lines are the fitted curves; (a) At/To = 
5.08 x 10 3 /(B tuh JB ) + 5.97 x 10, (b) At/r = 
-8.57 x (qH ) + 7.80 x 10, (c) max(V y )/C s0 = 
8.35 x 10" 3 x (B tuh JB ) - 3.26 x 10" 1 , (d) 
max(V y )/C s0 = S.10xlO- 1 x(qH )+S.02xlO- 2 . 



twist q is significantly small. The curvature is effective for the 
flux tube with the shortest wavelength A. Here we assume that 
the rising speed in the CZ is given as a terminal velocity where 
the buoyancy of the t ube equals the aerodynamic drag by the 
surrounding flow field (Par ker! 19751 iMoreno-Insertis & Emonetl 
119961) and the downward magnetic tension. Buoyancy, dynamic 
drag, and tension force acting on a unit cross-sectional area are 
written as 



/b = 



B 2 



SnHp' 



pVj_ 

^tube 



and 

h = 



B 2 
4ttR c 9 



(15) 



(16) 



(17) 



respectively, where H v = Hq(T/Tq) denotes the local pressure 
scale height, V z the tube's vertical speed, Cd the drag coefficient 
of order unity, and R c is the radius of curvature. The mechanical 
balance /b = /d + h yields the terminal velocity 



V z 



1 



1 



4C D p \2H V R c 



(18) 



First, let us discuss the curvature effect. In Equation (H"8lh 
the tension force is negligible for R c — > oo, while the ten- 
sion becomes effective when R c ~ 2H p . The relationship be- 
tween the curvature radius R c and the perturbation wavelength 
A is illustrated as Fig. [5^. Here, we write the tube's height as 
Az = Ztop(0 - ^tube- From this figure, we have 



R c - R c cos 6 = Az 
R c sin 6 = A 



which gives 

A= ^2R c Az-(Az) 2 . 



(19) 



(20) 



Thus, using the condition R c ~ 2H p , we obtain the critical wave- 
length for the tension to be effective: 



A c = ^AH v Az-{Az) 2 . 



(21) 



For instance, when the tube is halfway to the surface, i.e., 
Ztop/Ho = -50 and thus Az/Hq = 50, the local pressure scale 
height at this depth is H p /Hq ~ 21. Therefore, the critical wave- 
length is evaluated to be A c ~ 4l.2Ho, and the flux tube with a 
wavelength smaller than this value will be resisted by the ten- 
sion force, A ^ A c ~ 41.2Hq. In Fig.Q]:, we found that only the 
tube with A = 25Ho shows slower emergence due to the effective 
curvature force, which satisfies the condition A ^ A c . 

Next, let us go on to the dependencies on the field strength 
and the twist, by considering R c — > oo. Now the equation of the 
terminal velocity (TTST) reduces to 



^^^^beV^exp^ 
e~ 2 R tubQ 



a 

tube 



8C D /f p p 



(22) 



Here, in the first line we use Equation (IT2l) and in the second line 
we assume r ~ Rtub&- From this equation, we see that the rising 
velocity is in simple proportion to the initial field strength B tube 
when q is constant. If we change q with considering 7? tu b e = 5//o, 
for qHo = [0.05, 0.1, 0.2], the third term in the right-hand- side 
of Equation (l22t gives 



4 



l+^ube- 




(23) 



That is, the second term has only a weak positive correlation to 
the value of q. Therefore, the rising velocity of the flux tube is 
proportional to the field strength, while it is almost independent 
on the initial twist. The trend of the rising speed found in Section 
13. H is thus explained. 
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flux tube 




surface 




Fig. 5. (a) Cross-section of the rising flux tube along the axis (in 
the x - z plane), (b) Cross-section of the rising flux tube and the 
plasma layer ahead of the tube (in the y - z plane). 



It should be noted h ere that L and q are not always lin early 
correlated. According to Mor eno-Insertis & Emo net (119961) . the 
boundary of the expanded tube is well defined by the equipar- 
tition surface, where the kinetic energy density equals the mag- 
netic energy density of the azimuthal field: 



1 



Bi 



2 pV = 



(25) 



On the basis of an analogy from Equation (TT2K the profile of 
the expanded tube at the equipartition surface, where the radial 
distance is r\ (the subscript "1" indicates the expanded tube), is 
assumed to be written as 



B x (n) = £ t ubeiexp|- — 



B<p(n) = qnB x (n) 



tubel 



(26) 



where # t ubei and R tu bei are the axial field and the typical radius, 
respectively. Then, Equation (l25t reduces to 



n 



exp 



^tubel 



( 4npV 2 



tubel 



1/2 



1 



Rtubei 



(27) 



Here we assume that V is more or less approximate to V z and 
thus V does not depend on q. Other values of 2?tubei> ^tubei, and 
p are also assumed to be constant for different q. Here we in- 
troduce the notation \i = n/R^^i. Then Equation (1771) reduces 
to 



yuexp(-yu )=l/q, 



(28) 



Note here that, if we substitute H p ~ 40//o, p ~ 275po (val- 
ues at ztube), B tuhQ = 67 B , q = 0.1 /H , and R tuhQ = 5H (val- 
ues for case A) and assume Cd ~ 1 in Equation (l22k we ob- 
tain Voo = 0.21C s o, which is comparable to the simulation result 
~ 0.17C s o (Fig.[T]). This agreement indicates that Equation (l22l) 
is a rather reasonable estimation of the tube's risin g speed (see 
also iParkerll 19751 iMoreno-Insertis & Emonel 1996k . 



4.2. Dependence on the Twist 



We found in Section [331 that when the twist q is varied while the 
field strength # tu b e is kept constant, the duration of the HDF At is 
almost constant while the maximum horizontal speed max(V y ) 
is proportional to q. 

This feature can be explained by considering a simple model 
illustrated as Fig.0). Here the flux tube with a head size of 2L is 
rising at V z , which pushes the plasma layer with a thickness D. 
The thickness D is also described as D ~ |ztop(0l> where z to p(0 = 
-Ztube +J(^ V z (f)df. From the discussion in Section |4Jl V z and thus 
D are independent of q, which indicates that the HDF duration 
At - D/V z is also independent of q. 

If we write the outflow speed as V y , mass flux conservation 
can be written as 



y - Y±l 



(24) 



Here, V z /D is independent of q. The head size of the tube L, 
however, depends on the twist q, since the aerodynamic drag 
peels away the tube's outer flux and its amount depends on the 
twist. The head size remains larger with q, which results in the 
stronger HDF; V y oc L(q). Thus the maximum speed, max(y y ), 
will also depend on q. 



or, 

.2 



jj - In /j = In cj, 



(29) 



where q = Rtubei [^ ubel / (4npV 2 )] l/2 . For a larger radial dis- 
tance, r\ » Rxxfoti, i.e., \i » 1, 



li ~ In q. 
Then we obtain 



n ~ ^tubel 



1/2 



(30) 



(31) 



Therefore, the effective size of the expanded tube L ~ r\ is at 
least positively correlated to q, but not in a linear manner. 

4.3. Dependence on the Field Strength 

When the field strength at the tube's axis # t ube is varied while 
the twist q is fixed, the maximum HDF speed, max ( Vy), is found 
to be roughly proportional to the field strength. From Equation 
(l24lh if we assume L/D is constant, the horizontal speed V y and 
thus the maximum speed max (V y ) are proportional to the field 
strength £ tube . 

As for the HDF duration At, however, Fig.Hk clearly shows 
two regimes: stronger field cases that show a decreasing trend, 
and a middle case that deviates from this trend. Thus we should 
take into account the difference between these regimes. 

First, let us focus on the stronger field regime. Since the 
rising speed is proportional to the field strength, stronger tubes 
emerge faster. In this case, the accumulated plasma ahead of the 
tube does not drain down so much because of the short emer- 
gence period, and thus the thickness of the plasma layer becomes 
almost the same for these cases. That is, the thickness of the layer 



6 



S. Toriumi & T. Yokoyama: 3D MHD simulation on the solar flux emergence 



D is constant and is independent of the rising speed V z . Since the 
rising speed is proportional to the field strength # t ube, we have 

At - D/V z ocl/B tuhe . (32) 

Hence, the HDF duration is inversely proportional to the field 
strength, which explains the trend in the stronger field regime of 
the fitted inverse function in Fig. [4^. 

As for the middle strength case, the emergence takes longer 
and thus the drainage of the accumulated plasma becomes more 
effective, resulting in the much thinner layer D compared to the 
rising speed V z . Therefore, the HDF duration At D/V z be- 
comes shorter and thus deviates from the inverse trend in the 
stronger field cases. 
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5. Summary 

In this parametric survey, we vary the axial field strength, twist, 
and perturbation wavelength of the initial flux tube. As a result, 
we found the following features. 

The rising speed in the CZ strongly depends on the initial 
field strength but its correlation with the twist is weak. The 
emergence is resisted by the curvature force only in the case 
when the perturbation wavelength is shortest. According to the 
analytic study, the rising rate (terminal velocity) is written as 

Vco oc # tube ^1 + q 2 R^ ubQ , which indicates a strong dependence 

on the field strength and a weak correlation with the twist. 

As the flux tube approaches the surface, the accumulated 
plasma ahead of the tube escapes horizontally around the sur- 
face layer, which is called the HDF. The driver of the HDF is 
found to be the pressure gradient. 

When the field strength increases, the maximum HDF speed 
becomes higher, because the rising speed mainly depends on the 
field strength. The HDF duration is divided into two groups. For 
the stronger tube regime (# t ube ~ 1002?o), me duration is in sim- 
ple inverse proportion to the field strength, while the weaker field 
regime (# t ube ~ 100#o) deviates from the trend in the stronger 
tube regime because the fluid draining is more effective. 

The duration of the HDF is found to have no relation with the 
tube's twist, since the rising speed is independent of the twist. 
However, the maximum HDF speed shows a positive correlation 
with the twist. This feature is explained by considering the head 
size of the main tube that remains after the aerodynamic drag 
peels away the tube's outer field. The head size remains larger 
with the twist, which results in the stronger HDF. 

If we apply the above dependencies of the HDF to the ac- 
tual observations, we may be able to obtain information on the 
magnetic field in the subsurface layer, which we cannot observe 
optically. 
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